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Abstract. We apply the Smaller Alignment Index (SALI) method to a 4-dimensional mapping of accelerator dynamics in 
order to distinguish rapidly, reliably and accurately between ordered and chaotic motion. The main advantage of this index 
is that it tends exponentially to zero in the case of chaotic orbits, while it fluctuates around non-zero values in the case 
of quasiperiodic trajectories. Thus, it avoids the notorious ambiguities concerning the eventual convergence of (maximum) 
Lyapunov exponents to (positive) non-zero values. Exploiting the different behavior of SALI in these two cases we produce 
phase space 'charts' where regions of chaos and order are clearly identified. Evaluating the percentage of chaotic and escaping 
orbits as a function of the distance from the origin we are able to estimate rapidly and accurately the boundaries of the 
dynamical aperture of a proton beam, passing repeatedly through an array of magnetic focusing elements. 



INTRODUCTION 

One of the basic problems in accelerator physics is 
the determination of the stability regions (islands of 
quasiperiodic motion) of a particle beam, as well as its 
dynamical aperture, i.e the domain about the ideal circu- 
lar path in which the particles' motion remains bounded 
[1]. In the case of 'flat' hadron beams, where the horizon- 
tal (x-) motion is dominant, 2-dimensional (2D) area- 
preserving mappings are frequently used to model the 
effect of nonlinearities as the particles repeatedly pass 
through focusing magnetic elements [2], or experience 
the beam-beam interaction with particles of a colliding 
beam [3]. The presence of invariant curves around the 
origin of such mappings (corresponding to the particle's 
ideal circular orbit), guarantees the long-time stability 
of the beam. In that case, the chaotic motion that exists 
between these invariant curves remains always bounded 
and so the beam particles do not escape to infinity. On 
the other hand, in the case of 'elliptical' bunched hadron 
beams, where the vertical (y-) and longitudinal (z-) mo- 
tion cannot be neglected and which are modelled by 
k = 2M-dimensional symplectic mappings with M > 1, 
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chaotic regions can be connected, providing a path for 
the orbits to move away from the origin and eventually 
escape to infinity. This practically means the loss of par- 
ticles in the storage rings of the accelerator and hence a 
reduction of the beam's dynamical aperture. 

In the present paper, we face the problem of the es- 
timation of the dynamical aperture of a 4D symplectic 
mapping, which describes the motion of a hadron parti- 
cle as it passes through a magnetic focusing element of 
the FODO cell type. This system has already been stud- 
ied in [4, 5, 6, 7]. In particular Vrahatis et al. [7] tried 
to approximate invariant tori existing as far away from 
the origin as possible by computing stable periodic or- 
bits of very high period. Their idea was that the corre- 
sponding islands of stability and the quasiperiodic 'tori' 
around these periodic orbits could play the role of an ef- 
fective barrier of orbital diffusion, although their pres- 
ence does not exclude the possible 'leaking' of chaotic 
orbits to large distances from the origin. 

A more direct approach to the problem of estimating 
the size of the dynamical aperture, is the actual charac- 
terization of orbits on a fine grid around the origin as or- 
dered or chaotic (escaping or not escaping to infinity). In 
this way one can determine the region about the ideal cir- 
cular path where predominantly ordered orbits exist, thus 
guaranteeing the stability of the beam, at least up to the 
number of iterations that the orbits have been computed. 



This approach requires the use of a fast and reli- 
able method which can distinguish between ordered and 
chaotic motion rapidly, reliably and accurately. The usual 
method of the computation of the maximal Lyapunov ex- 
ponent [8, 9] does not meet these requirements as the 
number of iterations needed for the convergence of the 
Lyapunov exponent to its actual (zero or positive) value 
is not known a-priori and may become very high. There- 
fore, the application of this procedure to a huge num- 
ber of initial conditions becomes impractical and its ef- 
ficiency as a tool for studying the dynamical aperture of 
particle beams seems seriously doubtful . For these rea- 
sons, we prefer apply to our system the Smaller ALign- 
ment Index (SALI) method [10, 11, 12, 13], which has 
been proved to be an efficient, reliable and very fast 
method of chaos detection. 



THE SALI METHOD 

The SALI method was introduced in [10] and has already 
been applied successfully to distinguish between ordered 
and chaotic motion in various mappings and Hamiltonian 
systems [11, 12, 13], as well as problems of Celestial 
Mechanics [14, 15], Galactic Dynamics [16, 17], Field 
Theory [18] as well as non-linear lattices [19, 20]. 

In order to recall the definition of the SALI let us con- 
sider the fe-dimensional phase space of a conservative dy- 
namical system, e. g. a 2M-dimensional symplectic map- 
ping or a Hamiltonian flow of N degrees of freedom, with 
k = 2N. In a symplectic mapping the evolution of an or- 
bit with initial condition X(0) = {x\ (0),x 2 (0), . . .x k (0)), 
k = 2M is governed by the discrete-time equations of the 
mapping, having the form 

X(n + l)=F(X(n)), (1) 

where X(n) = (x\ (n),X2(n), . . . ,Xk{n)) is the orbit's lo- 
cation at the n-th iteration of the mapping. On the other 
hand, in a Hamiltonian flow, the motion of an orbit with 
initial condition X(0) = (x l (0),x 2 (0), . . .x k (0)), k = 2N 
is governed by Hamilton's equations of motion, which 
have the general form 

d -^=F(X(t)), (2) 

where X(t) — (x\ (t ),X2(t), ■ ■ ■ ,x k (t)) is the orbit's posi- 
tion in the phase space at time t . 

Suppose we wish to determine the chaotic vs. ordered 
nature of an orbit X{n) of a symplectic mapping (orX(f) 
of a Hamiltonian system) with initial condition X(0). 
To do so, one traditionally follows the evolution of one 
deviation vector V(n) = (dxi(n),dx2(n), . . .,dx k {n)) (or 
V(f)), which can be considered as initially pointing to 
an orbit nearby the one under study, and computes the 



orbit's maximal Lyapunov exponent [8, 9]. In the case of 
mappings the evolution of such a deviation vector V(n) 
is governed by the equations of motion of the so-called 
tangent map: 

V(n + l)=DF(X(n))-V(n), (3) 

while, in the case of Hamiltonian flows we use the set 
of linear ordinary differential equations called variational 
equations: 

^l=DF(X(t)).V(t), (4) 

where DF denotes the Jacobian matrix of equations (1) 
or (2) evaluated at the points of the orbit under study. 
Since, in the present paper, we study symplectic map- 
pings, our notation from now on will be restricted to dif- 
ference equations, although the following concepts can 
be easily extended to the case of differential equations 
describing Hamiltonian flows. 

For the evaluation of the SALI we follow the time evo- 
lution of two initially different deviation vectors Vi(0), 
V 2 (0) and define SALI [10] as: 

SALI(n)=min{||yi(n)+t/2(n)||,||t/i(n)-t/ 2 (n)||} (5) 

where || • || denotes the usual Euclidean norm and V,-, 
i = 1,2 are normalized vectors with norm equal to 1, 
i. e.V i (n)=V i (n)/\\V i (n)\\. 

Two different behaviors of SALI are thus distin- 
guished: 

1 . If the orbit under study is chaotic, the two vectors 
Vi(n), V 2 (n) tend to coincide (or become opposite) 
along the most unstable direction (corresponding 
to the maximal Lyapunov exponent). In this case, 
SALI tends exponentially to zero following a rate 
which depends on the difference between the two 
largest Lyapunov exponents [13]. 

2. If the orbit is ordered (quasiperiodic), there is no 
unstable direction and vectors Vi(ra), V 2 (rc) tend to 
become tangent to the corresponding torus, having 
in general different directions. In this case SALI re- 
mains different from zero, fluctuating around some 
mean value [12]. An exception to this behavior ap- 
pears for ordered orbits of 2D mappings where the 
SALI tends to zero following a power law (see [10] 
for more details). 

The simplicity of SALI's definition, its completely dif- 
ferent behavior for ordered and chaotic orbits and its 
rapid convergence to zero in the case of chaotic motion 
are the main advantages that make SALI an ideal chaos 
detection tool, perfectly suited for multidimensional con- 
servative systems, such as proton (antiproton) beams in 
accelerator storage rings. 



GLOBAL DYNAMICS OF A 4D 
ACCELERATOR MAPPING 



Consider the 4D symplectic mapping 
/ 
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describing the instantaneous sextupole nonlinearities ex- 
perienced by the dynamics of a proton beam as it passes 
repeatedly through magnetic focusing elements of the 
FODO cell type [4, 5, 6, 7]. Here x\ and X3 are the par- 
ticle's horizontal and vertical deflections from the ideal 
circular orbit of the beam, x 2 and X4 are the associated 
momenta and 01 and (02 are related to the accelerator's 
tunes q x , q y by 



a>i = 2nq x , ©2 = 2nq y 



(7) 



Let us first examine the behavior of the S ALI for some 
individual orbits. Vrahatis et al. [7] have computed near 
the boundary of escape of the mapping several stable pe- 
riodic orbits of very long period, as well as some invari- 
ant tori near them. In Figure 1(a) we see the projection 
into the x\ —X2—X3 space of one such ordered orbit of 
the mapping (6) for q x = 0.61803 and q y = 0.4152 first 
presented in figure 1(b) of [7]. This orbit is generated by 
a small perturbation in q x of the stable periodic orbit of 
period 13237 found in [7] for q x = 0.61903 and it lies on 
8 tori in the x\ — xi — xj, space. The exact values of the 
orbit's initial condition, which will be denoted from now 
on as X*(0) = (x*(0),4(0),x?j(0),4(0)), can be found 
in table 3 of [7]. Following [8, 9], we may compute the 
maximal Lyapunov exponent X\ , of the orbit as the limit 
for n — > 00 of the quantity 



Li(n) 
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n ln ||V(0)|| 



i. e. Ai = lim L\ (n), (8) 



where V (0), V (n) are deviation vectors from the orbit at 
n = and n > iterations respectively. We recall that 
A; = for ordered orbits, while Ai > for chaotic orbits. 
The ordered nature of the orbit is clearly revealed by the 
evolution of its SALI (Figure 1(b)) and of its L\ (Figure 
1(c)). The SALI remains different from zero fluctuating 
around 10~ 2 5 converging to this value long before L\ 
becomes zero. Thus, as far as SALI is concerned, the 
computation could have been stopped after about 10000 
iterations, concluding correctly that the orbit is ordered. 




(6) 





FIGURE 1. An ordered orbit of mapping (6). (a) Projection 
of 50000 consequents of the orbit into the x\ — xi — xj space. 
Evolution of orbit's SALI (b) and L\ (c) as a function of 
mapping's iterations n in log-log scale. 



By changing the x\ coordinate of the initial condition 
of the orbit of Figure 1, to x\ — 0.624 we get the weakly 
chaotic orbit plotted in Figure 2(a), surrounding the 8 tori 
of Figure 1(a). This behavior is clearly seen in Figure 
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FIGURE 2. A chaotic orbit of mapping (6). (a) Projection 
of 50000 consequents of the orbit into the x\ — X2 ~ x 3 space, 
(b) Projection on the x\ — X2 plane of the points of panel (a) 
(black points) and of Figure 1(a) (gray points) with |*3 1 < 0.04. 
Evolution of orbit's SALI (c) and L\ (d) as a function of 
mapping's iterations n in log-log scale. 



2(b) where we project on the x\ — X2 plane the points 
of Figures 1(a) and 2(a) with |x 3 < 0.04. The SALI of 
the chaotic orbit decreases rapidly to zero, after a tran- 
sient initial phase (Figure 2(c)) reaching the limits of the 
computer's accuracy (i. e. 10~ 16 ) after about n = 20000 
iterations, showing clearly the chaotic nature of the or- 
bit. We can of course set a less demanding threshold for 
the SALI's values in order to define an orbit as chaotic. 
Considering for example as such a threshold the value 
SALI w 10~ 8 , we can characterize the orbit as chaotic 
after only 8000 iterations. For the same number of iter- 
ations, L\ (Figure 2(d)) does not seem to converge to a 
non-zero value, so that many more iterations are needed 
for the definite characterization of the orbit as chaotic by 
the use of the maximal Lyapunov exponent. In fact, at 
about 10000 iterations the maximal Lyapunov exponent 
gives an erroneous picture, as it starts to fall to values 
closer to zero! 

In Figure 3(a) we see the projection into the x\ — 
X2 — X3 space of a chaotic orbit near the edge of the 
beam's dynamical aperture, with initial condition X*(0) 
forq x = 0.628615, q y = 0.4152, which escapes to infinity 
after about n = 82000 iterations. Again, the SALI rapidly 
determines the chaotic nature of the orbit as it becomes 
less than 10~ 8 after about n = 12000 iterations (Figure 
3(b)), while L\ continues to decrease showing no sign 
of convergence to a non-zero value, until after 32000 
iterations (Figure 3(c)). 

This fundamentally different behavior of the SALI for 
ordered (Figure 1(b)) and chaotic orbits (Figures 2(c) and 
3(b)) and its rapid determination allows us to perform 
efficiently a more 'global' study of the dynamics of 
mapping (6) in order to estimate the region of stability 
around the origin. 

As a first step in that direction, let us compute, up 
to n = 10 s iterations, a great number of orbits whose 
xi(0) coordinate varies from to 0.9 with a step equal 
to 10~ 4 , while X2(0), *3(0), *4(0) are the same as in the 
stable orbit of Figure 1 . In Figure 4 we plot the SALI 
of each orbit (after n = 10 s iterations) as a function of 
the initial coordinate *i(0). We note that chaotic orbits 
which escape in less than n — 10 s iterations are charac- 
terized as escaping orbits and are not plotted in Figure 
4. From Figure 4 we see that orbits with x\ (0) < 0.54 
are ordered, having SALI > 10~ 4 , except for a tiny inter- 
val around x\ (0) « 0.39 where one finds chaotic orbits 
having SALI w 10~ 9 . The region 0.54 < x x (0) < 0.59 is 
occupied mainly by chaotic orbits having SALI values 
less than 10~ 8 and is followed by a region of mainly or- 
dered motion for 0.59 < x\ (0) < 0.79. For larger values 
of x\ (0) chaos dominates while escaping orbits that are 
not plotted in Figure 4 appear for x\ (0) > 0.83. 

Since we are interested after all in global picture of the 
dynamics around the origin, it is convenient to use the 
SALI method for 'charting' this region. Let us consider 
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FIGURE 3. A chaotic escaping orbit of mapping (6). (a) 
Projection of 82000 consequents of the orbit into the x\ — X2 — 
*3 space. Evolution of orbit's SALI (b) and L\ (c) as a function 
of mapping's iterations n in log-log scale. 

therefore orbits with initial conditions on a grid mesh 
around the origin and evolve them for a given number n 
of iterations. We shall characterize each orbit as chaotic 
if SALI < 1(T 8 and as ordered if SALI > lfT 8 . If the 



FIGURE 4. The values of the SALI for n = 10 5 iterations of 
orbits with constant initial coordinates X2(0), x?,(0), xn(0) and 
x\ (0) e [0, 0.9] , as a function of x\ (0) . 



orbit escapes before the final number n of iterations is 
reached it will be characterized as an escaping orbit. 

We first restrict our study to the 2-dimensional con- 
figuration plane x\ — xj, in order to be able to visual- 
ize our results. In particular, we consider orbits on a 
400 x 400 grid mesh uniformly covering the rectangu- 
lar region —1 < xi(0) < 1, —1 < ^3(0) < 1, keeping 
fixed thex 2 (0), X/t(0) values. The corresponding 'charts' 
are plotted in Figure 5 for X2(0) = xzt(0) = (Figures 
5(a),(b)) and forx 2 (0) =x 4 (0) = 0.1 (Figures 5(c),(d)). 
The orbits were followed for n = 10 4 (Figures 5(a), (c)) 
and n = 10 5 (Figures 5(b), (d)) iterations in order to un- 
derstand the evolution of 'charts' with respect to the 
number of iterations. In Figure 5 the initial conditions 
corresponding to chaotic orbits are plotted in black while 
the ones corresponding to ordered orbits are plotted in 
gray. In all panels of Figure 5 all non-colored points cor- 
respond to escaping orbits. 

From the comparison of panels (a) and (b), and pan- 
els (c) and (d) of Figure 5 we see that the region occu- 
pied by non-escaping orbits (ordered and chaotic ones) 
does not practically change as the number n of iterations 
increases. This means that most of the escaping orbits 
fly away from the central region very fast, after a small 
number of iterations. So, the initial conditions plotted by 
black and gray color in Figure 5 define the region of sta- 
bility around the beam's circular motion in the sense that 
all these orbits do not escape. We also see that in all pan- 
els of Figure 5 the region around the origin corresponds 
to ordered motion, while chaotic orbits exist mainly at 
the borders of the stability region. As the number of it- 
erations increases, the number of chaotic orbits also in- 
creases. This happens because weakly chaotic orbits lo- 
cated at the borders of the region of ordered motion re- 
veal their chaoticity later on as their SALI needs more 
iterations in order to become less than 10~ 8 . Thus, al- 
though the number of non-escaping orbits remain practi- 
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FIGURE 6. Percentages of ordered, escaping (solid curves) 
and chaotic (dashed curves) orbits, within spherical shells of 
width dr = 0.01 as a function of shell's mean radius, r ce n, 
from the origin. The characterization of the orbits as chaotic 
or ordered was done according to their SALI values after n = 
10 4 (gray curves) and n = 10 5 (black curves) iterations. The 
percentages of escaping orbits do not change significantly as 
the number of iterations increases and so the corresponding 
curves practically coincide. 
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FIGURE 5. Regions of different values of the SALI on the 
x\ — X3 plane after n = 10 4 iterations (panels (a) and (c)) 
and after n = 10 iterations (panels (b) and (d)). The initial 
conditions of the computed orbits on the X2 — X4 plane are 
^2(0) = ^4(0) = for panels (a) and (b) and ^2(0) = X4(0) = 
0.1 for panels (c) and (d). In all frames, initial conditions 
are colored black if their SALI < 10~ 8 and gray if SALI > 
10~ 8 . The uncolored initial conditions correspond to orbits that 
escape in less than n iterations. 



cally constant the percentage of this number that corre- 
sponds to chaotic orbits increases as n grows. 

Considering orbits with initial conditions uniformly 
distributed around the origin within a 'volume' of the 
full 4-dimensional phase space, we now perform a more 
global analysis of orbital stability. As we cannot produce 
plots like the ones of Figure 5 for the 4-dimensional 
space, we present in Figure 6 the percentages of the 
various types of orbits within spherical shells of width 
dr = 0.01 inside a 4-dimensional hypersphere of radius 
r = 1 centered at the origin. We note that by the distance 
r of an initial condition (x\ (0),X2(0),X3(0),X4(0)) from 
the origin (0,0,0,0) we refer to the quantity 

r = ^/4(0)+x2(0)+x2(0)+x2(0) _ (9) 

From the results of Figure 6 we see again that the number 
of escaping orbits does not change significantly as the 
number n of iterations increases, while the percentage 
of chaotic orbits increases with n. An estimation of the 
radius of the dynamical aperture therefore gives r w 0.55, 
as up to that distance from the origin no escaping orbits 
are found. Of course, for 0.4 < r < 0.55 there exists a 
significant amount of non-escaping chaotic orbits. So 
r f» 0.4 is a reasonable estimate of the maximal radius 
of a 4-dimensional hypersphere around the origin, where 
orbits not only do not escape to infinity but, in addition, 
are also ordered. 



CONCLUSIONS 

In the present paper, we have applied the method of 
the Smaller Alignment Index (SALI) to the character- 
ization of orbits of a 4D symplectic mapping describ- 
ing the dynamics of a proton beam passing repeatedly 
through magnetic focusing elements of the FODO cell 
type. Thus, we have been able to locate efficiently is- 
lands of ordered motion, layers of weak chaos, as well as 
estimate accurately the beam's dynamical aperture. 

The success of this approach lies in the fact that it 
can rapidly distinguish between ordered and chaotic mo- 
tion in Hamiltonian flows and symplectic mappings of 
any dimensionality. Since the SALI decays exponen- 
tially to zero in the case of chaotic orbits (and oscillates 
quasiperiodically around non-zero values in ordered re- 
gions), it bypasses the slow and often irregular conver- 
gence properties of the computation of Lyapunov expo- 
nents and thus provides quickly a definite characteriza- 
tion of each orbit. 

This allows one to use the SALI to characterize whole 
domains in phase space of different scales and sizes and 
'chart' them as chaotic or regular. Carrying out such a 
study for the mapping of this paper, we have been able 
to 'trace out' the dynamical aperture of proton beams 
with a 2-dimensional (x and y) cross section, by locat- 
ing 4-dimensional domains, where non-escaping behav- 
ior is guaranteed even after a very high number of it- 
erations. Currently, we are extending our work to more 
realistic 6-dimensional mappings, where longitudinal (or 
synchrotron) oscillations are taken into consideration and 
space charge effects are included [21]. Despite the addi- 
tional complications present in these models, we believe 
that the SALI method will again be able to yield useful 
results, 'charting' correctly the dynamics of phase space 
domains that would otherwise be very difficult to probe 
efficiently by more traditional techniques. 
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